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The recently proposed diagrammatic expansion (DE) technique for the full Gutzwiller wave func¬ 
tion (GWF) is applied to the Anderson lattice model. This approach allows for a systematic 
evaluation of the expectation values with full Gutzwiller wave function in the finite dimensional 
systems. It introduces results extending in an essential manner those obtained by means of stan¬ 
dard Gutzwiller approximation (GA) scheme which is variationally exact only in inhnite dimensions. 
Within the DE-GWF approach we discuss principal paramagnetic properties and their relevance to 
the heavy fermion systems. We demonstrate the formation of an effective, narrow /-band originat¬ 
ing from atomic /-electron states and subsequently interpret this behavior as a direct itineracy of 
/-electrons; it represents a combined effect of both the hybridization and the correlations reduced 
by the Goulomb repulsive interaction. Such feature is absent on the level of GA which is equivalent 
to the zeroth order of our expansion. Formation of the hybridization- and electron-concentration- 
dependent narrow /-band rationalizes common assumption of such dispersion of / levels in the 
phenomenological modeling of the band structure of GeGoIns. Moreover, it is shown that the 
emerging /-electron direct itineracy leads in a natural manner to three physically distinct regimes 
within a single model, that are frequently discussed for 4/- or 5/- electron compounds as separate 
model situations. We identify these regimes as: (i) mixed-valence regime, (ii) Kondo-insulator bor¬ 
der regime, and (iii) Kondo-lattice limit when the /-electron occupancy is very close to the /-states 
half-filling, (n/) ^ 1. The nonstandard features of emerging correlated quantum liquid state are 
stressed. 

PACS numbers: 71.27.-|-a, 71.10.-w, 71.28.-|-d, 71.10.Fd 


I. INTRODUCTION AND MOTIVATION 

Heavy fermion systems (HFS) belong to the class of 
quantum materials with strongly correlated 4/ or 5/ 
electrons. They exhibit unique properties resulting from 
their universal electronic features (e.g. very high den¬ 
sity of states at the Fermi level) almost independent 
of their crystal structure. Among those unique prop¬ 
erties are: (i) enormous effective masses in the Fermi- 
liquid state, as demonstrated through the linear specific 
heat coefficient^— and their direct spin-dependence in 
the de Haas-van Alphen measurements^^— , (ii) Kondo- 
type screening of localized or almost localized /-electron 
magnetic moments by the conduction electron o^d^ , (iii) 
unconventional superconductivity, appearing frequently 
at the border or coexisting with magnetisroii, and (iv) 
abundance of quantum critical points and associated with 
them non-Fermi (non-Landau) liquid behavior— 

The Anderson lattice model (ALM), also frequently 
referred to as periodic Anderson model, and its deriva¬ 
tives: the Kondoi^Ai and the Anderson-Kondo ^^d^ lat¬ 
tice models, capture the essential physics of HFS. Al¬ 
though, the class of exact solutions is known for this 
modelSSTi^^, they are restricted in the parameter space. 
Thus, for thorough investigation of the model properties 
the approximate methods are needed. One of the earli¬ 
est theoretical approaches for the models with a strong 
Coulomb repulsion was the variational Gutzwiller wave 


function (GWF) method^i^— . However, despite its sim¬ 
ple and physically transparent form, a direct analytic 
evaluation of the expectation values with full GWF can¬ 
not be carried out rigorously for arbitrary dimension and 
spatially unbound systems. 

One of the ways of overcoming this difficulty is the 
so-called Gutzwiller Approximation (GA), in which only 
local two-particle correlations are taken into account 
when evaluating the expectation values. GA provides 
already a substantial insight into the overall properties 
of strongly correlated systems^ d°d^d°~ — . Moreover, this 
approach has been reformulated recently to the so-called 
statistically-consistent Gutzwiller approximation (SGA) 
scheme and successfully applied to a number of problems 
involving correlated electron system a^^d^~ — . Among 
those, a concrete application has been a microscopic de¬ 
scription of the fairly complete magnetic phase diagram 
of UGe 2 ^i^ which provided quantitatively correct re¬ 
sults, even without taking into account the 5/-orbital 
degeneracy due to uranium atoms. 

An advanced method of evaluating the expectation val¬ 
ues for GWF is the variational Monte Garlo technique 
(VMG)ii-.^. However, this method is computationally 
expensive and suffers from the system-size limitations. 
Though, one must note that the VMG method allows 
for extension of GWF by including e.g. Jastrow intersite 
factors^. 

Here we use an alternative method of evaluating the 
expectation values for GWF, namely a systematic dia- 
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FIG. 1: (Color online) Upper part: Schematic representation 
of the two-orbital Anderson lattice model with initially local¬ 
ized /- and delocalized c-electrons, and hybridization between 
them. Bottom part: Emerging quasiparticle states in the hy¬ 
bridized bands of correlated particles. On the left: the shapes 
of the density of states in the respective situations. 


grammatic expansion for the Gutzwiller wave function 
(DE-GWF)^“— . This method was formulated initially 
for the Hubbard model in two dimensions in the context 
of Pomeranchuk instability^, and applied subsequently 
to the description of high-temperature superconductivity 
for the Hubbard^i^ and the t-J— models. In the zeroth 
order of the expansion this approach straightforwardly 
reduces to the GA^. For the one-dimensional Hubbard 
model it converges^ to the exact GWF results. Within 
DE-GWF a larger variational space can be sampled than 
within the alternative VMC technique because the long- 
range components of the effective Hamiltonian are ac¬ 
counted for naturally. The DE-GWF method (truncated 
to match the variational spaces) reproduces the results 
of VMC with improved accuracy (as shown for the t-J— 
and the Hubbard models^). Additionally, the method 
works also in the thermodynamic limit. In effect, the ap¬ 
proach is well suited to capture subtle effects, e.g. those 
related to the topology of the Fermi surface in the cor¬ 
related stated or the investigated here formation of a 
narrow /-electron band. 

In this study, we extend the DE-GWF approach to dis¬ 
cuss principal paramagnetic properties within ALM. The 
emergence of the quasiparticle picture is schematically il¬ 
lustrated in FiglU Explicitly, we investigate the shape of 
the quasiparticle density of states (DOS, p{E)), evolving 
with the increasing order of the expansion, k. For fc > 0 
the hybridization gap widens up with respect to that in 


GA {k = 0 case) and DOS peaks are significantly pro¬ 
nounced. Moreover, we investigate the DOS at the Fermi 
level {p{Ep)) evolution with the increasing the hybridiza¬ 
tion strength \V\ ~ total electron concentration n, plane, 
as it is a direct measure of the heavy-quasiparticle effec¬ 
tive mass. We find that this parameter is significantly 
enhanced for fc > 0, mainly in the low hybridization 
limit and at the border of the Kondo-insulating state. 
Furthermore, we trace the contribution coming from the 
originally localized /-electrons (cf. Fig. [T]- upper part) 
to the quasiparticle spectrum with the increasing order 
of the expansion. For fc > 0, /-quasiparticles effectively 
acquire a nonzero bandwidth (up to 6% of the conduction 
bandwidth) as a combined effect of both interelectronic 
correlations and hybridization. 

Assumption of a narrow / band existence has recently 
been made in a phenomenological modeling of the heavy 
fermion compound CeGoIns band structure^^— . We 
show that the emergence of such a band, absent in GA 
(fc = 0), is an evidence of the /-electron direct itineracy 
explained later. To quantify this itineracy we introduce 
the parameter Wf - the width of the effective, narrow /- 
band. On the hybridization strength - total electron con¬ 
centration, |V| - n plane, Wf is significantly enlarged in 
the three distinct regimes, which we identify respectively 
as the mixed-valence, Kondo/almost Kondo insulating, 
and the Kondo-lattice regimes (when /-electron concen¬ 
tration is close to the half-filling, i.e., when (fif) —>■ I). 
These physically distinct regimes are frequently discussed 
and identified in various experiments^ d^’^^~ — and in 
theoryiSi^^iSliSS. 

The structure of the paper is as follows. In Sec. El we 
describe the ALM Hamiltonian and define the Gutzwiller 
variational wave function in a nonstandard manner. In 
Sec. mil we derive the DE-GWF method for ALM and 
determine the effective single-particle two-band Hamilto¬ 
nian. In Sec. m we present results concerning param¬ 
agnetic properties: the quasiparticle spectrum, the resul¬ 
tant density of states at the Fermi level, and formation 
of an effective narrow /-electron band out of initially lo¬ 
calized states. In Appendix 1X1 we discuss the equivalence 
of the zeroth-order DE-GWF approach with GA. In Ap¬ 
pendix E] we present some technical details of DE-GWF 
technique. 


II. MODEL HAMILTONIAN AND 
GUTZWILLER WAVE FUNCTION 

Our starting point is the Anderson lattice model 
(ALM) with the chemical potential p and expressed 
through Hamiltonian 

- Y + Y^^f~ 

iJ.<T i.<T i.a 

EUY^f^L + + VJclM, 

i iJ,o' 
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where i = (ix, iy) (and similarly j) is the two-dimensional 
site index, (//o-) (^L) annihilation 

(creation) operators related to /- and c- orbitals respec¬ 
tively, and a =t, i is the z-component direction of the 
spin. We assume that the hopping in the conduction 
band takes place only between the nearest neighboring 
sites, tij = the hybridization has the simplest 

onsite character—, Vij = the local Coulomb re¬ 

pulsion on the / orbital has the amplitude U, and the 
initially atomic / states are located at the energy e/. In 
the following |t| is used as the energy unit. 

Gutzwiller wave function (GWF) is constructed from 
the uncorrelated Slater determinant |'0o) by projecting 
out fraction of the local double /-occupancies by means 
of the Gutzwiller projection operator Pq, 

\^g) = Pg\^o) = X{Pg-M- ( 2 ) 

i 

In the GA approach when only a single / orbital (in the 
present case) is correlated the projection operator can be 
defined by 

Pg\ = 1 - (1 - g)ni^n{^, ( 3 ) 

where g is a variational parameter. Such form allows 
for interpolating between the fully correlated {g = 0) 
and the uncorrelated {g = 1) limits. Equivalently one 
can consider average number of doubly occupied states, 
= d? as a, variational parameter. 

The Gutzwiller projection operator can be selected dif¬ 
ferently as proposed in Ref. [TOj, namely 

pl,pG-i=p^,i = i+xdr. (4) 

In the above relation a; is a variational parameter and for 
the paramagnetic and translationally invariant system we 
define Hartree-Fock (HF) operators of the form 

= (nft - nof){h(^ - no/), (5) 

where n^f denotes average occupation of a single / 
state and spin a in the uncorrelated state, |'0o), be., 
’^o/ = {fiafiPo- Hereafter the shortened notation for 
the expectation values is used, i.e., ('!/'o|---|'0o) = 
Strictly speaking, although, d^^ has not the Hartree- 
Fock form of the double occupancy operator, the HF su¬ 
perscript has its meaning as the property, {df^)o = 0 is 
preserved. 

On the other hand, the Gutzwiller projection operator 
can be defined in general form as 

= |r)i(r|i, (6) 

r 

with variational parameters Ar G {Aq, A/, A/, A^} 
that characterize the possible occupation probabil¬ 
ities for the four possible atomic Fock /-states 


Relation (|4]) couples Ar and x, reducing the number 
of independent variational parameters to one. Explicitly, 
we may express the parameters Ar by the coefficient x, 

A^ = l + xnlf, 

Pi=Pi = >i = ^-xnof{l-nof), ( 7 ) 

A^ = 1 -|- x(l — riQ/)^- 

As the parameters Ar and x are coupled by the condi¬ 
tions 0, there is a freedom of choice of the variational 
parameter; in this work we have selected x. The pa¬ 
rameter X covers the same variational space as g in GA. 
Additionally, the projector (g]) leads to much faster con¬ 
vergence than ([3]) (cf. Ref. [s^. From dH) it is clear 
that X = 0 corresponds to the uncorrelated limit. The 
other extremity, the fully correlated state is reached for 
X = max{a:(Ad = 0),a;(Ao = 0)}. This leads to the 
bounds max{} < x < 0. The minimal 
value is X = —4 for ng/ = 0.5. 

The method is suitable for an arbitrary filling of the 
/ orbital. However, due to the fact that present work is 
mainly addressed to the description of the Ce-based com¬ 
pounds, we study the regime in which the /-orbital fill¬ 
ing either does not exceed unity or is only slightly larger. 
Precisely, in the all figures presented here the /-orbital 
filling is never larger than 1.05. 

III. DE-GWF METHOD 
A. General scheme 

In this section we present general implementation of 
the DE-GWF method. The procedure is composed of 
the following steps: 

1. Choice of initial state \'ipo). 

2. Evaluation of {U)g = selected \tpo) - 

cf. Sec, mm 

3. Minimization of {'H)g with respect to the varia¬ 
tional parameter (here x). 

4. Construction of the effective single particle Hamil¬ 
tonian determined by 6pl^^{\tpo)) = ^?7(|^/’o)) - cf. 

Sec. rrrrci 

5. Determination of \iPq) as a ground state of the ef¬ 
fective Hamiltonian - cf. Sec. IlHDl 

6. Execution of the self-consistent loop: starting again 
from the step 1 with |'0 q) until a satisfactory con¬ 
vergence, i.e., \iPq) = IV'o)) is reached. 

Steps 4 and 5 ensure that the final form of |'0o) rep¬ 
resents the optimal choice which minimizes the ground 
state energy {H)g- The DE-GWF method with respect 
to other related methods, GA and VMG, introduces a 











4 


new technique for evaluating the expectation value of the 
correlated Hamiltonian with GWF (step 2 of the above 
procedure). In particular, it provides an important im¬ 
provement as, e.g., for GA only single sites in the lattice 
contain the projection whereas the remaining environ¬ 
ment does not. GA leads e.g. to the inability of obtaining 
the superconducting phase in the Hubbard model^. On 
the other hand, the VMC method tackles that problem 
properly, but at the price of extremely large computing 
power needed. This leads to the lattice size limitations 
(typically up to 20x20 sites) and a limited distance of 
real space intersite correlations taken into account. 

In this respect, DE-GWF introduces, in successive or¬ 
ders of the expansion, correlations to the environment of 
individual sites (beyond GA), as well as converges in a 
systematic manner to the full GWF solution. Also, DE- 
GWF was shown to provide results of better accuracy 
than VMC—, and additionally, is free from the finite- 
size limitations. It also demands definitely less computa¬ 
tional power than VMC. Thus in general, this method is 
capable of treating more complex problems with GWF. 
On the other hand, DE-GWF is tailored specifically for 
GWF, while VMC allows for starting from different forms 
of variational wave function e.g., adding the Jastrow 
factors^i^. 


B. Diagrammatic expansion 

The key point of the variational procedure is the cal¬ 
culation of the expectation value of Hamiltonian © with 
GWF I^g) (point 1 from the scheme in Sec. IHI Al) . by 
starting from the expression 

^ ii’G I Ti I Ipc) ^ ii’O I Pg'HPg I '0o) .O'! 

^ - (V'G I ^g) {^0 \P^\^0) ' ^ ^ 

We use the DE-GWF technique^^— , based on the 
expansion of the expectation values appearing in Eq. 

in the power series in variational parameter x, 
with the highest power representing number of corre¬ 
lated vertices assumed to be correlated in the envi¬ 
ronment - besides local ones. This method is sys¬ 
tematic in the sense that the zeroth order corresponds 
to GA— , whereas with the increasing order the full 
GWF solution is approached. Explicitly, we deter¬ 
mine expectation values with respect to GWF of any 
product operator originating from the starting Hamil¬ 
tonian (P) C>i(j) = 

is executed by first accounting for the projection 
part on the site i(j) - external vertices (e.g., computing 
= pGAPG;j)Oi(^j){pG-j)pG-i, see below) and then, 
including one-by-one correlations (terms) to the other 
sites 1 7^ i, j - internal vertices. 

Formally, the procedure starts in effective power ex¬ 


pansion in X of all relevant expectation values 

(v^G I Oio) I v-g) = n Ph\ 

oo k (9) 

k=0 

where The prime in the multiple 

summation denotes restrictions: Ip ^ Ip/, and Ip 7^ i,j for 
all p,p'. k is the order of the expansion. Note that for 
k = 0 we obtain {tjjc\Oi(j)\tjjG) = This means 

that the projection operators act only locally (i.e., only 
the sites i and j are affected) and in this case we recover 
the GA results (for a details discussion of the equiva¬ 
lence see Appendix Expectation values in can 
now be calculated by means of the Wick’s theorem in 
its real-space version, as they involve only products av¬ 
eraged with \ijjo)- Such power expansion in x allows for 
taking into account long-range correlations between k in¬ 
ternal sites (li, ..., Ifc) and the external ones (i, j). It must 
be noted that it is not a perturbative expansion with re¬ 
spect to the small parameter x. Instead, the expansion 
should be understood as an analytic series with the order 
determined by the number of correlated internal vertices 
taken in the nonlocal environment. For k = 00 , the full 
GWF solution would be obtained. However, on the basis 
of our results, a satisfactory results for the expansion in 
ALM case are reached already starting from fc = 3. 

As said above, the expectation values (...)o in Eq. m 
can be evaluated by means of the Wick’s theorem. Then, 
the terms with k internal sites can be visualized as dia¬ 
grams with k internal and 1 (or 2) external vertices. The 
lines connecting those vertices are defined as, 

P'U = (Cio-Cj(T)oi 

= (^lcj.)o, (10) 

Tij = {fiafiPo — Sijnof. 

By constructing the projector operator we have 
eliminated all the diagrams with the local /-orbital 
contractions {{fi^fi<T)o), the so-called Hartree bubbles. 
This procedure, as discussed explicitly in Ref. [H, leads 
to significantly faster convergence than that with the 
usual Gutzwiller projector, with the variational param¬ 
eter g^. It constitutes the main reason for the effi¬ 
ciency of the DE-GWF method. Finally, all the ex¬ 
pectation values with respect to GWF are normal¬ 
ized by (^/gIV’g) (cf. Eq. ®). However, through the 
linked-cluster theorem^^, the terms coming from ex¬ 
pansion of ('0g|'0g) = ('0o|.PgI^o) cancel out with all 
disconnected diagrams appearing in the numerator of 
Eq. ([8]). In effect, the expectation values can be ex¬ 
pressed in the closed form by the diagrammatic sums 
S G jc(2) j/(2) j/(4)|^ defined 

in Appendix IbJ what leads to the following resultant ex- 
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pression for the ground state energy: 


D. Determination of |i/)o) 


Wg 

L 




+2(e/ - /i) ^no/ + (1 + 

+C/A^ (do + 2nofI^^^^ + (1 - 

+W + /3T^'=(3,i)^ ^ 


( 11 ) 


where the trivial sums ~ Si = -^ have already 

been included. Parameters {a,/3, 7 , m, do} are all func¬ 
tions of no/and x (cf. AnnendixlBl Eg. (IB2I) '). For fc = 0, 
only the diagrammatic sums /c( 2 ) 2 ^/o(i,i) 

not vanish and we reproduce the standard GA result; the 
Coulomb energy reduces to UX^do and hybridization to 
4Pa(/lci)o, whereas the diagrammatic sums for c-band 
only are trivial (cf. the discussion in Appendix [AT) . 

The expectation value {'H)g calculated diagrammat- 
ically is minimized next with respct to the variational 
parameter x (step 3 in the scheme in Sec. IIII Al) . 


C. Effective quasiparticle Hamiltonian 


In this section we determine \tpQ) as a ground state of 
(point 5 from the scheme in Sec. IIII Ap . 

In order to obtain the effective dispersion relations for 
c- and /-electrons and the k-dependent hybridization we 
use the lattice Fourier transform 


Af) _ 1 V pdi-j)k/c(/) 
- L 2^^ bj ’ 


IJ 


= i V pdi-j)kv/c 

k t/-^ ij ■ 


(15) 


1 .J 


In this manner, we reduce the many-body problem to 
the effective single-quasiparticle picture (cf. Fig. [T|) de¬ 
scribed by the effective two-band Hamiltonian. The 2x2 
-matrix representation of Eq. (US resulting from such a 
transform, has the following form 









E("L fL)r^ 

\r rr ^ 



(16) 


The next step in our procedure (step 4 in the scheme in 
Sec. IIII Al) is the mapping of the correlations contained 
in {iPg\T~L\'4!g) ! {'4>g\'‘Pg) onto the corresponding uncorre¬ 
lated expectation value (/’o|’H®®|'!/’o)- It is realized via 
the condition that the minima of the expectation values 
of both Hamiltonians coincide for the same equilibrium 
values of lines (nni) and no/, which define |/’o)- Note that 
the present formulation of this step of our minimization 
procedure is equivalent to those previously used^“—. Ex¬ 
plicitly, 


W®)o(C', F, IT, no/) = (5 (H)g(C', E, IT, no/) 


d{n)c 

dC 


■5C- 


d{n)c 

dW 


AIT- 


d{n)c 

OF 


■SF- 


d{H)c 

dnof 


Ano/, 


( 12 ) 

where skipping lattice indices for lines means that we 
consider each of them separately. It leads directly to the 
following form of the effective single-particle two-band 
Hamiltonian with non-local interband hybridization, i.e, 


+ Y t(JL& 

ij,(7 

ij,<7 


(13) 


where the effective hopping and hybridization parameters 
are derivatives with respect to lines. 


A- — 

ly - 


d{n)c 

dCu 


d{n)G 

OF, ’ 


y_fc ^ d{n)G 

« aiTij ^ 
d{'H)G 


(14) 


tl = 


dn 


Of 


where the eigenvalues, of the above Hamiltonian are 

I^k. = C + «v/(C? + ^^, (17) 

where a = ±1 differentiates between the two hybridized 
bands. For convenience, we have defined 

e+ = and (18) 

F in Eq. (11611 is the unitary transformation matrix to 
the basis in which is diagonal, defined as 



It is now straightforward to obtain the principal correla¬ 
tion functions (lines), i.e. 

(4aCkfT)o = ■u+0(Fk+) -k Mi0(Fk_), 

(/k<TCkfT)o = u+u- (0(Fk+) - 0(Fk-)), (21) 

(/L/k-)o = n 2 _ 0 (Fk+) +n^ 0 (Fk-), 

where 0(F) denotes the Heaviside step function and 
plays the role of an energy cutoff for respective quasipar¬ 
ticle bands energies (II3- Using the reverse Fourier trans¬ 
formation we obtain self-consistent equations for lines 
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and no/, 

ka 

W"iJ = ;^E(/k.5k.)o eW 

( 22 ) 

^iJ = lE(/kJk.)o eW, 

ka 

’^o/ = -j- E^(/ka-/kg-)o- 

ka 

To determine the properties of the model, we solve 
in the self-consistent loop the system of Eqs. ca and 
(point 6 from the scheme in Sec. IIII All 

Finally, the ground state energy Eq is defined by 

Eg = ('H)gIo + (23) 

where (H)g|o denotes the expectation value (fTTll of the 
starting Hamiltonian for the equilibrium values of the 
lines and the total number of particles is defined by 
n = 2{n(^ + nf^)G- The /-orbital filling separately is de¬ 
fined by n/ = 2(h(jG- 

IV. RESULTS AND DISCUSSION 

A. System description and technical remarks 

In our analysis we consider a square, translationally 
invariant, and inhnite {L —>■ oo) lattice, with two orbitals 
(/ and c) per site. The square lattice consideration is 
justified by the common quasi-two-dimensional layered 
structure of / atoms in the elementary cell of many Ce- 
based heavy fermion systems^di that our studies are rel¬ 
evant to. 

While proceeding with the diagrammatic expansion 
(DE), in principle two approximations need to be made. 



10 12 3 

IJx-'xl 

FIG. 2: Schematic illustration of the real-space cutoff on the 
lattice. The solid lines denote exemplary, in terms of dis¬ 
tance, correlation functions (referred to as lines) taken into 
account between i-site (in the center) and the j-sites (on the 
periphery). Farther connections are not considered. 


First, only the lines (flO)) satisfying the relation |i — jp = 
iix—jx)'^ + {,iy—jy)'^ < 10 are taken into account (i.e., we 
make a real-space cutoff - cf. Fig. [ 21 ). For comparison, 
in VMC only rarely lines farther than these connecting 
nearest neighboring sites (more precisely, only the lines 
corresponding to the hopping term range of the starting 
Hamiltonian) are taken into accoun t . From our nu¬ 
merical calculations it follows that the nearest- and the 
second-nearest neighbor contractions compose the domi¬ 
nant contributions (cf. Fig. | 71 d ). 

The second limitation in DE is the highest order of the 
expansion, k, taken into account. Asymptotic behavior 
starting from fc = 3, of some properties such as the den¬ 
sity of states (DOS) at the Fermi level (FL), p{Ef), and 
the width of the effective / band, Wf (cf. FigsjH jS}: and 
speak in favor of the calculation reliability, achieved 
already in that order. Therefore, if not specified oth¬ 
erwise, the expansion is carried out up to the third or¬ 
der ik = 3), i.e., with the three internal vertices taken 
into account. We stress again that the zeroth-order ap¬ 
proximation (fc = 0) is equivalent to the GA approach 
(cf. Appendix [^ for details). The results of GA are 
regarded here as a reference point for determining a sys¬ 
tematic evolution, including both qualitative and quan¬ 
titative changes, when the higher-order contributions are 
implemented. 

The parameters of the ALM Hamiltonian m are taken 
in units of |t|: a strong Goulomb repulsion is taken as 
U = 10, the reference energy for /-electrons, e/ = —3, 
the onsite hybridization is assumed negative and varies 
in the range |U| G (0.8, 2.5), and the total band filling 
(n = 2{n(^ + nf^)G) is in the range allowed by the condi¬ 
tion that the / level occupancy per site {rif = 2{n(^)G) 
roughly does not exceed unity. The reason for considera¬ 
tion of this regime is the circumstance that for interesting 
us Ge-based compounds the concentration of / electrons 
per cerium should not exceed 1 (i.e., with the Ge^’*' and 
Ge^"*" configurations only). However, from the construc¬ 
tion of the method the regime for n/ > 1 is fully accessi¬ 
ble and physically correct. In carrying out the DE-GWF 
procedure we adjust the chemical potential p = Ep for 
the fixed total filling n. Numerical integration of Eq. 
(1^^ and the self-consistent loop were both performed 
with precision of the order of 10“® or better with the 
help of Gnu Scientific Library (GSL) procedures'^. 


B. Correlated Fermi liquid 

Before the detailed analysis is carried out, a method¬ 
ological remark is in place. The effective Hamilto¬ 
nian (fl^ is of single-particle form, but coupled to the 
self-consistent procedure of evaluating the relevant av¬ 
erages (1^ . However, this does not compose the full 
picture. The physical quantities are those obtained 
with a projected wave function. For example, nj = 
E^(V'G|/iL/i<T|^/’G) = Y.a{EGn{^'PG)o: which in general 
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FIG. 3: Difference between uncorrelated and correlated /- 
electron number, An/ = along the line 

of constant hybridization, |I/| = 1.3, with respect to chang¬ 
ing total filling. The specific character of the region IV is 
explained in Sec. EYl 

is slightly different from X]ct('™L)o- The situation is il¬ 
lustrated explicitly in Fig. [S] In effect, the quasipar¬ 
ticle picture is amended with the nonstandard features 
of this correlated (quantum) liquid (CL). Parenthetically, 
the same difference will appear when considering mag¬ 
netic and superconducting states, where the magnetic 
moments, {S?)a vs. {S?)o, and the superconducting 
gaps, (Aij)G and (Aij)o will be different. So, we have a 
mapping of the correlated onto quasiparticle states, but 
not of the physical properties. In brief, we have to dis¬ 
tinguish between the correlated and the uncorrelated /- 
electron occupancy or other property even though, from 
the way of constructing CSD, the density of quasiparticle 
states (coming from (TT^ b represents that in the corre¬ 
lated state. 


C. Quasiparticle Density of States 

We start with analysis of the quasiparticle DOS emerg¬ 
ing from the DE-GWF method in successive orders of the 
expansion (cf. Fig. 0 ]). For k > 0 and the total filling 
n = 1.97 (i.e., near the half filling), the hybridization 
peaks become more pronounced (cf. Fig0]-the inset Ta¬ 
ble) and the hybridization gap increases. 

For k > 0 the overall shape of DOS changes only quan¬ 
titatively (cf. Fig. 0]). However, the value of the DOS at 
the Fermi level, p{Ep), changes remarkably (cf. the inset 
to Fig. 0|). Although for fc = 1 it is underestimated and 
for k = 2 overestimated, for fe = 4 we see no significant 
difference with respect to the fc = 3 case. For this rea¬ 
son, if not specified explicitly, the subsequent analysis is 
proceeded in the third order, fc = 3. 

The value of p(Ef) is of crucial importance. This 
parameter is a measure of the quasiparticle effective 
mass, as the latter is inversely proportional to the second 



E-Ep 

FIG. 4: (Golor online) Density of states (DOS) near the half- 
filling (n = 1.97) at \V\ = 1.5 for selected orders of the di¬ 
agrammatic expansion (k — 0,3). Explicit values of p{Ef) 
are also listed in the inset Table (for 0 < fc < 4). For fc = 3 
a satisfactory convergence of the expansion is reached. The 
k = 1,2,4 plots are not included for clarity as, apart from 
peak heights, they are practically the same as the plot for 
fc = 3. For k > 0 (beyond GA) the hybridization peaks are 
more pronounced (large DOS at the Fermi level p(Ef)), which 
is related directly to the increased by correlation effective- 
mass enhancement for quasiparticles. 


derivative of the energy, at the Fermi surface, and 

thus is determined by p(Ep). 

In Fig. EK we draw the value of p(Ep) on the plane 
hybridization - total electron number (per site), V - n. 
This quantity is particularly strongly enhanced near the 
half filling (n ~ 2). In comparison to the lowest value 
p(Ep) Ri 0.75, the maximal enhancement is of the or¬ 
der of 40. In Fig. Eb we present evolution of p(Ep) 
on the logarithmic scale with the decreasing total filling 
and approaching n = 2 (vertical arrow in Fig. Et- marked 
by the encircled letter b ). The extrapolated value of 
p(Ep) may reach extremely high values of 1000 and even 
more (dashed line in Fig. Eb) in the region IV. Such fea¬ 
ture could explain extremely high mass renormalization 
in some of HFS for large but finite value of the Coulomb 
interaction U. 

The region where p(Ep) is enhanced strongly, is that 
with low hybridization \V\ values and for the total filling 
n ~ 1. This region is strictly correlated with the position 
of the second pronounced peak in DOS (cf. Fig. 0]) which 
therefore has its meaning as the Van Hove singularity. 
Additionally, for n/ ~ 1, where the effects of correlations 
are the strongest, we observe also a large value of p(Ep). 
In that limit the stability of magnetic phases should be 

studied separatel y 

As marked in Fig. Eb, near the total half-filling, n ~ 2, 
we could not obtain a satisfactory convergence of our 
self-consistent procedure. This is attributed to the po¬ 
sition of the chemical potential extremely close to the 
hybridization-induced peaks (significant when n/ > 0.9). 
Technically, this leads to extreme fluctuations (out of our 
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FIG. 6: (Color online) /-electron density of states Pf{E) 
within successive orders of expansion (fc < 4). For k = I 
and higher, formation of the effective /-band can be clearly 
observed. For fc = 3 the final shape of Pf{E) and the value of 
/-band width Wf stabilize. 


FIG. 5: (Color online) (a) Density of states at the Fermi level 
p{Ef) on the hybridization strength - total electron concen¬ 
tration plane, \V\ - n. Additionally (not marked), for n = 2 
we obtain always the Kondo insulating state. By IV (for con¬ 
sistency with FigO we have marked a V-shaped region where 
we have no numerical convergence due to the presence of sin¬ 
gular hybridization peaks for low | V| and with n near the half 
filling (see main text), (b) Evolution of p{Ef) in the half- 
logarithmic scale near the region IV (along the vertical arrow 
with the letter b). By extrapolation (dashed line in (b)), for 
the almost half-filled situation, p{Ef) can be enhanced even 
by factor of 1000 relative to its lowest values on the |V| - n 
plane, (c) Evolution of p{Ef) with the decreasing |V| (along 
the horizontal arrow with letter c), within successive orders 
of the expansion {k < 4). For large \V\ > 1.8, GA (fc = 0 
order) provides already realistic values of p{Ef)- 


numerical precision) of the filling, effective hopping pa¬ 
rameters, and the lines coming from the effective Hamil¬ 
tonian USD, as they are sensitive to a slight change of the 
chemical potential position. For n = 2 and nonzero hy¬ 
bridization, we obtain always the Kondo insulating state. 
However, strictly speaking, the true Kondo-type compen¬ 
sated state is demonstrated explicitly only if magnetic 
structure is taken into account explicitly^ dOd^ . 

In Fig. [St we depict the p{Ep) evolution with the de¬ 
creasing hybridization amplitude \V\ for fc < 4. Our 
results show that for large |K|, GA (fc = 0) already is 
reasonable approximation. The situation changes as we 
approach the low-|K| regime near the half-filling, where 
inclusion of higher-order contributions leads to a strong 
enhancement of p{Ef), as discussed above. 

In summary, the quasiparticle mass is enhanced spec¬ 
tacularly near n = 2 and in the regime of small hybridiza¬ 
tion \V\. The /-state occupancy is then n/ > 0.9. This 
is the regime associated with the heavy-fermion and the 
Kondo-insulating states. We discuss those states in detail 


in what follows. 


D. /-electron direct itineracy 

As stated already, the DE-GWF method is used here 
to map the correlated (many-body) system, described 
by the original Hamiltonian o with the help of the 
Gutzwiller wave function onto that described by 

the effective quasiparticle Hamiltonian (11311 with an un¬ 
correlated wave function IV'o)- By constructing the ef¬ 
fective Hamiltonian it is possible to extract the explicit 
contribution to the quasiparticle picture as coming from 
a direct hopping between the neighboring / sites. By 
contrast, in GA (fc = 0) case, the / electrons itineracy 
is only due to the admixture of c-states when the quasi¬ 
particle states are formed. Once we proceed with the 
diagrammatic expansion to higher order (fc > 0), they 
start contributing to the quasiparticle spectrum in the 
form of a dispersive /-band (cf. Fig. jH]). The resulting 
band is narrow, Wf < 0.5, whereas the starting conduc¬ 
tion (c) band has the width of Wc = 8. As was mentioned 
in the Sec. HI we interpret the parameter Wf as a measure 
of emerging degree of direct itineracy^ i.e., presence of a 
direct hoppings between the neighboring / states in the 
effective Hamiltonian. 

Again, a methodological remark is in place here on the 
numerical convergence of the results with respect to fc. 
Namely, the /-bandwidth appears already for fc = 1, but 
both its width and the curvature stabilizes only starting 
from fc = 3. 

In the recent phenomenological modeling of 
GeGoIng^^— the band structure used is the hybridized- 
two-band independent-particle model with dispersive 
/-band, even though the Ge 4/ states can be placed well 
above the so-called Hill limit, where there should not be 
any direct hopping between the original neighboring / 
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FIG. 7: (Color online) (a) Effective bandwidth of /-states, 
Wf, on the hybridization strength \V\ - electron concentra¬ 
tion n plane. Wf is regarded as a measure of direct itineracy 
of /-electron states. Three separate disjoint regions (light 
color) are regarded as universal and frequently discussed as 
separate limits, both in theory and experiment. Namely, 
the mixed valence regime (III), the Kondo/almost Kondo- 
insulating regime (II), and the Kondo-lattice regime (I) with 
ny ^ 1 — (5, 5 <C 1. (b) Effective /-electron intersite hop¬ 
pings ty along the marked vertical line of the diagram for 
\V\ = 1.3. The energy dispersion for / quasiparticles is deter¬ 
mined mainly by the nearest and the second nearest hoppings 
Region IV, near n = 2 is marked separately due to the 
lack of convergence of the numerical results (see main text). 


states. The fit presented there provides Wf of the same 
order of magnitude as that obtained here. As those 
phenomenological models do not include the Coloumb 
interaction, the ground state is determined by the 
uncorrelated wave function. Hence, our analysis of the 
effective Hamiltonian resulting from ALM provides a di¬ 
rect microscopic rationalization of the narrow dispersive 
f-band presence assumed ad-hoc in the fitting procedure 
in Ref. [H-Iwl. 

In Fig. [71i we display diagram comprising the width 
of /-band Wf on |H| ~ n plane, with contours of con¬ 
stant values of n/. We observe the appearance of re¬ 
gions, where the / quasiparticles have a sizable band¬ 
width (bright color) and other, where they remain local¬ 


ized (dark regimes). We expect that in the regions, where 
/ electrons are forming a band, a nontrivial unconven¬ 
tional superconductivity and/or magnetism may appear. 
These topics should be treated separately as they require 
a substantial extension of the present approach (incorpo¬ 
rating new type of lines 

With the help of the width Wf we may single out 
three physically distinct regimes (cf. Fig. [7^). We 
identify those regions as the mixed-valence regime (HI), 
the Kondo/almost Kondo insulating regime (H), and the 
Kondo-lattice regime (I) with Uf —>■ 1 —<5, with i5 <C I (cf. 
Fig. EJi). These universal regions are usually discussed 
independently within different specific models and meth¬ 
ods. In regime I the role of f-c Coulomb interactions 
(the Falicov-Kimball term) may be needed for complete¬ 
ness (cf. Ref. lli, whereas in the Kondo-lattice regime 
the transformation to the Anderson-Kondo model is ap¬ 
propriate (cf. Refs. [isllT^ . In the extreme situation, 
the heavy-fermion states are modeled by pure Kondo- 
lattice model2^“— . However, strictly speaking, the last 
model applies only in the limit of localized / electrons 
(rif = I), since then the total numbers of / and c elec¬ 
trons are conserved separately. 

In Fig. [71 d we present the effective hopping parame¬ 
ters for / states for \V\ = 1.3, i.e., along the marked 
vertical line in Fig. [7^. This line crosses three singled 
out regions of the itineracy. The leading contribution 
to the /-electron band energy arises from the nearest- 
and the second nearest-neighbor hoppings. Such circum¬ 
stance confirms that our earlier assumption about the 
real-space cutoff shown in Fig. [7] has been selected prop¬ 
erly. Moreover, it points to the importance of including 
also the components beyond those of the starting Hamil¬ 
tonian, only rarely taken into account within the VMC 
method^i^. 

In Fig. [5] we show the contributions to the effective 
hybridization. The initial (bare) local hybridization ac¬ 
quires momentum dependence. Nevertheless, the local 
part is still dominant since the nonlocal terms are at least 
two orders of magnitude smaller. 

The emerging in our model /-band introduces a new 
definition of the /-electron itineracy as it is not so much 
connected to the Fermi-surface size2^, but with the ap¬ 
pearance of a direct hoppings between / sites. This differ¬ 
ence is highly nontrivial, especially in the limit n/ = 1, 
where we obtain the largest bandwidth Wf. Such be¬ 
havior is attributed to the specific character of our ap¬ 
proach. Namely, we consider here the processes within 
our initial Hamiltonian o, but under the assumption 
that the neighboring sites are also correlated. This, as 
we have shown directly, leads also to the finite /-band 
in the effective single particle Hamiltonian (fT51) . The re¬ 
sults thus throw a new light on the longstanding issue 
of the dual localized-itinerant nature of / electrons in 
HFS ^9 i 80 . While the magnetism can be attributed to the 
almost localized nature of / electrons, an unconventional 
superconductivity requires their itineracy in an explicit 
manner, as will be discussed elsewhere^ 
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FIG. 8: (Color online) Spatial contributions to the effec¬ 
tive hybridization normalized by the first onsite (i = j) Vbo 
term along the marked vertical line of the diagram in FiglTja, 
for \V\ = 1.3. Note that due to correlation initially local, on¬ 
site hybridization acquires effectively momentum dependence. 
However, the nonlocal contributions constitute only up to 2% 
of the local one. 

V. SUMMARY 

We have applied a recently developed diagrammatic 
technique (DE-GWF) of evaluating the expectation val¬ 
ues with the full Gutzwiller wave function for the case 
of two-dimensional Anderson lattice. We have analyzed 
properties of the model by discussing the most important 
features of the heavy fermion systems in the paramagne¬ 
tic state. We have also shown that by approaching in suc¬ 
cessive orders of the expansion the full Gutzwiller-wave- 
function solution, we obtain a systematic convergence. 
In the zeroth order of expansion our method reduces to 
the standard Gutzwiller Approximation (GA). 

In difference with GA, DE-GWF does not overesti¬ 
mate the hybridization narrowing factor. Furthermore, 
our method produces unusually enhanced peaks at the 
Fermi level in the density of states, particularly near the 
half-filling, n —t 2. This in turn, is connected to the value 
of effective mass and by analyzing in detail this region we 
can explain a very large mass enhancement observed in 
heavy fermion systems as described by ALM with large, 
but finite Goulomb-interaction value, here U = 10|t|. 
The regions of sizable p{Ef) enhancement are also found 
in the small-hybridization limit and are connected to the 
presence of both the Van Hove singularity and the strong 
correlations in the limit of n/ —t 1. 

The /-electron contribution to the full quasiparticle 
spectrum is analyzed in detail. For nonzero order of the 
expansion {k > 0) we observe a systematic formation of 
the effective /-band with the increasing k. In spite of 
the fact that the bare electrons are initially localized, / 
quasiparticles contribute to the total density of states as 
they become itinerant. We interpret this property as the 
emerging direct /-electron itineracy. As a measure of this 
behavior, we introduce the the width Wf oi effective /- 


band. Formation of such narrow /-band rationalizes e.g. 
the recent phenomenological modeling of the CeGoIns 
band structure^^— . 

The nonstandard character of the resultant correlated 
Fermi liquid (GL) which differs from either the Landau 
Fermi liquid (FL) and the spin liquid (SL), should be 
stressed. FL represents a weakly correlated state (no lo¬ 
calization) and SL represents a fully correlated state. Our 
GL state in this respect has an intermediate character. 
Namely, the quasiparticle states are formed (as exempli¬ 
fied by e.g. density of states), but the physical properties 
such as the occupancy nj, the magnetic moment (S^) or 
the pairing gap in real space (Ajj) are strongly renormal¬ 
ized by the correlations. Such situation is often termed as 
that of an almost localized Fermi-liquid state^i^ d^d^d^ . 

By analyzing the results on the hybridization strength 
|V| - total band filling n plane, we single out explicitly 
three physically distinct regions, which we regard as three 
separate universality limits. Namely, we have linked 
those disjoint regions with the regimes frequently dis¬ 
cussed as separate classes in the heavy fermion systems: 
the mixed-valence regime, the Kondo/almost Kondo in¬ 
sulating regime, and the Kondo-lattice regime for n/ —>• 
1. We suggest, that the regions of significant /-electron 
itineracy can be connected to the unconventional heavy 
fermion superconductivity which would require separate 
studies. 

We have also commented on the longstanding issue of a 
dual localized-itinerant nature of / electrons in the heavy 
fermion systems. The new definition of itineracy is in 
accord with their (almost) localized nature. 
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Appendix A: Equivalence of the k=0 order DE-GWF 
expansion and the Gutzwiller approximation (GA) 

Here we show the equivalence of the zeroth order DE- 
GWF and the standard Gutzwiller approximation (GA). 
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In both methods (DE-GWF in the zeroth order of expan¬ 
sion k = 0) the effect of the projection can be summarized 
by the expressions for evaluating following expectation 
values: (hi'|-hi 4 ,)G and (//^Cio- -|- H.c.)g- The remaining 
averages in ALM are unchanged under the projection. 

Explicitly, in the DE-GWF for /c = 0 the resulting 
averages are expressed as follows 

(nitni4,)[?"°^ = Xlnlf (Ala) 

{fLci. + H.c.)[?=°) = aifZci. + H.c.)o, (Alb) 

where parameter a (see also Appendix [B] Eqs. (IBlI) and 
(IB2II 1 is defined as 

a = (1 — no/)AoAs-I-uo/AdAs- (A2) 

On the other hand, in GA the resulting averages are 
expressed as^ 


(A3a) 

{fZci, + H.c.)[f^) = + H.c.)o, (A3b) 

where the parameter (P is the double occupancy prob¬ 
ability, and q is the so-called Gutzwiller factor reducing 
the hybridization amplitude, which for the equal number 
of particles for each spin is defined as 


^ _ Z ~ ~ ^^0/ + d‘^) + V (^0/ ~ 

\/no/(l - no/) 


(A4) 

If we identify double occupancy probabilities expressed 
by both methods in (IA3a|) and (lAlal) to be equal, yielding 
cP = A^ngy, then the parameter a (jA2p exactly reduces 
to the parameter ^ (IA4I) . 

GA procedure results in the effective single-particle 
Hamiltonian of the form 


'^GA = 

k,cr 


V^v\ 


^k<T + LUd'^ 


-^n ( Y “ ^” 0 /) -^^rn(Y ' 

k.tr k,cr 

(A5) 

In the above Hamiltonian it is necessary to add con¬ 
straints for /-electron concentration and their magnetiza¬ 
tion in order to satisfy consistency of the procedure^LS^. 
In effect, the whole variational problem is reduced to 
minimization of the ground state energy with respect to 
d^, no/, mf, and respective Lagrange multipliers A)) and 
A^, playing the role of the effective molecular fields^. 
However, the effect of constraint for /-electron magne¬ 
tization is relevant only in the case of magnetism con¬ 
sideration either as intrinsi o'^^d^ or induced by applied 
magnetic field^. Here, as we discuss paramagnetic state 

= >'L = Hi = 0- 

The DE-GWF method by construction guarantees that 
the variationally obtained /-electron occupancy number 
n/ coincides with that obtained self-consistentlyS. We 


k=0 k=l 


k=2 


-pCCd.l) 

' ij 
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FIG. 9: (Color online) Diagrammatic sums to the second or¬ 
der, k — 2. c- and /- orbital sites are denoted with empty and 
filled circles respectively. Solid, dashed (blue) and dotted (or¬ 
ange) connections represent F, W, and C lines respectively 
(cf. Eq. (1101) 1. The numbers in brackets under diagrams 
stand for their multiplicity resulting from the Wick’s theo¬ 
rem. Note that by construction of our sums we have no di¬ 
agrams with so-called “Hartree bubbles”, namely loop-lines 
within the same site and orbital. 


have thus provided analytical argument for the equiva¬ 
lence of the DE-GWF method for fc = 0 and the standard 
GA procedure. Also, by an independent numerical cross¬ 
check we have verified that all the observables calculated 
within both methods indeed coincide. 


Appendix B: Diagrammatic sums 


We start with expressions for the following projected op¬ 
erators originating from ALM Hamiltonian ([TJ, namely 


PoZiPc-,! = Ad[2no/hf^ + (1 - xdo)dP^ + doP^i], 
PG-,inicrPo-i = (1 -h xm)nP^ -h ydf ^ -h nofPY 

pG-,jHpG:i = afH + Pflaur, 

(Bl) 
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where additionally we have defined 


with the fc-th order contributions 


nHF ^ f^HF ^ f^HF 

/3 = Xs{Xd — Aq), 

OL = AfiAo + f3nof, 
7 = a:(l - 2nof), 


do — ^ 0 /) 

m = no/(l - no/). 


(B2) 


As mentioned in the main text, such form of the pro¬ 
jected operators significantly speeds up the convergence 
of the numerical results^, since by construction all two- 
operator averages for a single site and /-orbital, the so- 
called Hartree bubbles, vanish. The above operator alge¬ 
bra leads to the compact definition of the diagrammatic 
sums: S-e 
Eq. (fTTI) . 


OO ^ 

^ = T.^S{k). (B3) 


J,cc(l,l)(^) 

i^‘^\k) 

i^^\k) 


(4aCjadq^...i,)S, 


E os. 


/o- 




(B4) 


Superscript c in the expectation values means that only 
the connected diagrams are to be included. Note that 
in (IB4I) there are no summation restrictions, due to the 
linked cluster theoremi^. The resulting diagrammatic 
sums for S up to second order (k = 2) are depicted in 
Fig. O 
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